1 Outline





2 Results


set.seed(12345) # for reproducibility
options(knitr.kable.NA = '')

# install // load packages 
# Some packages need to be loaded. 
# We use `pacman` as a package manager, which takes care of the other packages. 
if (!require("distill", quietly = TRUE)) install.packages("distill")
if (!require("devtools", quietly = TRUE)) install.packages("devtools")
if (!require("papaja", quietly = TRUE)) devtools::install_github("crsh/papaja")
if (!require("patchwork", quietly = TRUE)) devtools::install_github("thomasp85/patchwork")
if (!require("klippy", quietly = TRUE)) devtools::install_github("RLesur/klippy")
if (!require("pacman", quietly = TRUE)) install.packages("pacman")
if (!require("Rmisc", quietly = TRUE)) install.packages("Rmisc")
if (!require("rstatix", quietly = TRUE)) install.packages("rstatix")
if (!require("effsize", quietly = TRUE)) install.packages("effsize")
if (!require("lsr", quietly = TRUE)) install.packages("lsr")
if (!require("effectsize", quietly = TRUE)) install.packages("effectsize")
if (!require("ggbeeswarm", quietly = TRUE)) install.packages("ggbeeswarm") # Never load it directly.
pacman::p_load(tidyverse, papaja, knitr, dplyr, car, psych, afex, lme4, lmerTest, 
               emmeans, ggplot2, ggpubr, lattice, latticeExtra, parallel, effects, psycho, apa)
library("patchwork"); library("klippy")
klippy::klippy()


2.1 Phase 1 - Categorization


t1 <- read.csv("data/reactMEM_p1_all.csv", header = T)
glimpse(t1, width = 70)
## Rows: 32,640
## Columns: 15
## $ SN      <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ Trial   <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 1…
## $ Block   <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ bTrial  <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 1…
## $ cCong   <int> 1, 1, 1, 1, 2, 2, 2, 2, 1, 1, 2, 2, 2, 2, 1, 1, 2, 2…
## $ cnRep   <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2…
## $ cnPair  <int> 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2…
## $ cIdt    <int> 1, 1, 2, 2, 3, 3, 4, 4, 5, 5, 6, 6, 4, 4, 5, 5, 3, 3…
## $ cCat    <int> 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, 1, 2, 2, 1, 2, 1, 2, 1…
## $ Im_nExm <int> 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1…
## $ Im_Idx  <int> 199, 193, 47, 77, 157, 127, 209, 111, 181, 72, 61, 1…
## $ Im_Name <chr> "199_3_out.jpg", "193_1_out.jpg", "047_3_in.jpg", "0…
## $ Resp    <int> 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, 1, 2, 2, 1, 2, 1, 2, 1…
## $ RT      <dbl> 1.1000, 1.0478, 0.7386, 0.7283, 0.8903, 0.9315, 0.70…
## $ Corr    <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
length(unique(t1$SN))
## [1] 32

# change class of main factors: double to factor
t1$SN = factor(t1$SN)
t1$Block = factor(t1$Block)
t1$cCong = factor(t1$cCong, levels=c(1,2), labels=c("v","nv")) 
t1$cnRep = factor(t1$cnRep, levels=c(1,2,3,5,4), labels=c(1,2,3,4,0))
t1$cnPair = factor(t1$cnPair)
t1$cIdt = factor(t1$cIdt)
t1$cCat = factor(t1$cCat)
t1$Im_Name = factor(t1$Im_Name)
t1$Im_Idx = factor(t1$Im_Idx)
t1$RT <- t1$RT*1000;
t1$Corr <- as.numeric(t1$Corr==1)

# check number of trials for each condition/SN
table(t1$SN)
## 
##    1    2    3    4    5    6    7    8    9   10   11   12   13   14   15   16 
## 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 
##   17   18   19   20   21   22   23   24   25   26   27   28   29   30   31   32 
## 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020 1020
table(t1$Block, t1$SN) 
##    
##       1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19
##   1 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##   2 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##   3 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##   4 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##   5 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##   6 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170 170
##    
##      20  21  22  23  24  25  26  27  28  29  30  31  32
##   1 170 170 170 170 170 170 170 170 170 170 170 170 170
##   2 170 170 170 170 170 170 170 170 170 170 170 170 170
##   3 170 170 170 170 170 170 170 170 170 170 170 170 170
##   4 170 170 170 170 170 170 170 170 170 170 170 170 170
##   5 170 170 170 170 170 170 170 170 170 170 170 170 170
##   6 170 170 170 170 170 170 170 170 170 170 170 170 170
table(t1$cCong, t1$SN) 
##     
##        1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18
##   v  540 540 540 540 540 540 540 540 540 540 540 540 540 540 540 540 540 540
##   nv 480 480 480 480 480 480 480 480 480 480 480 480 480 480 480 480 480 480
##     
##       19  20  21  22  23  24  25  26  27  28  29  30  31  32
##   v  540 540 540 540 540 540 540 540 540 540 540 540 540 540
##   nv 480 480 480 480 480 480 480 480 480 480 480 480 480 480
headTail(t1)
##         SN Trial Block bTrial cCong cnRep cnPair cIdt cCat Im_nExm Im_Idx
## 1        1     1     1      1     v     1      1    1    2       3    199
## 2        1     2     1      2     v     1      2    1    2       1    193
## 3        1     3     1      3     v     1      1    2    1       3     47
## 4        1     4     1      4     v     1      2    2    1       1     77
## ...   <NA>   ...  <NA>    ...  <NA>  <NA>   <NA> <NA> <NA>     ...   <NA>
## 32637   32  1017     6    167    nv     0      2  118    1       3     54
## 32638   32  1018     6    168     v     4      2  119    2       1    144
## 32639   32  1019     6    169    nv     4      2  120    1       1     13
## 32640   32  1020     6    170    nv     4      2  118    1       1     54
##             Im_Name Resp     RT Corr
## 1     199_3_out.jpg    2   1100    1
## 2     193_1_out.jpg    2 1047.8    1
## 3      047_3_in.jpg    1  738.6    1
## 4      077_1_in.jpg    1  728.3    1
## ...            <NA>  ...    ...  ...
## 32637  054_3_in.jpg    1 1392.9    1
## 32638 144_1_out.jpg    2  765.9    1
## 32639  013_1_in.jpg    1 1465.3    1
## 32640  054_1_in.jpg    1 1380.5    1
glimpse(t1, width = 70)
## Rows: 32,640
## Columns: 15
## $ SN      <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ Trial   <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 1…
## $ Block   <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ bTrial  <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 1…
## $ cCong   <fct> v, v, v, v, nv, nv, nv, nv, v, v, nv, nv, nv, nv, v,…
## $ cnRep   <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 2…
## $ cnPair  <fct> 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2, 1, 2…
## $ cIdt    <fct> 1, 1, 2, 2, 3, 3, 4, 4, 5, 5, 6, 6, 4, 4, 5, 5, 3, 3…
## $ cCat    <fct> 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, 1, 2, 2, 1, 2, 1, 2, 1…
## $ Im_nExm <int> 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1, 3, 1…
## $ Im_Idx  <fct> 199, 193, 47, 77, 157, 127, 209, 111, 181, 72, 61, 1…
## $ Im_Name <fct> 199_3_out.jpg, 193_1_out.jpg, 047_3_in.jpg, 077_1_in…
## $ Resp    <int> 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, 1, 2, 2, 1, 2, 1, 2, 1…
## $ RT      <dbl> 1100.0, 1047.8, 738.6, 728.3, 890.3, 931.5, 701.8, 8…
## $ Corr    <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…


2.1.1 RT analysis: General



## 1.2.1. preprocessing
# filter null/incorrect data
ct1 <- t1 %>% filter(Corr == 1) # filter(Corr == 1) # remove incorrect trial

# check accuracy
100-100*(nrow(ct1)/nrow(t1))
## [1] 4.666054
# 5.76% incorrect trials were not analyzed

# trimming 3sd outlier trials
# tt1 <- ct1
tt1 <- ct1 %>% filter(RT > 200 & RT < 10000) %>%
  group_by(SN) %>% # grouping by participants
  nest() %>%
  mutate(lbound = map(data, ~mean(.$RT)-3*sd(.$RT)),
         ubound = map(data, ~mean(.$RT)+3*sd(.$RT))) %>% # make new data (3sd cut)
  unnest(c(lbound, ubound))%>%
  unnest(data) %>%
  mutate(Outlier = (RT < lbound)|(RT > ubound)) %>% # set outlier
  filter(Outlier == FALSE) %>% # filtering outlier
  ungroup() %>%
  dplyr::select(SN, Block, cCong, cnRep, cnPair, cIdt, cCat, Im_Name, RT, Corr) # select variables to analyze

# outlier trial ratio
100-100*(nrow(tt1)/nrow(ct1))
## [1] 1.812514
## [1] 1.89857

# RTs were trimmed
# when 1) faster than 200ms 2) slower than 10 secs 3) 3SD away from the participants mean

# mean number of trials for each conditions
tt1 %>% group_by(SN, cCong, cnRep) %>% 
  filter(cnPair == 2) %>% 
  summarise(NumTrial = length(RT)) %>%
  ungroup() %>%
  group_by(cCong, cnRep) %>%
  summarise(Mean = mean(NumTrial), 
            Median = median(NumTrial), 
            Min = min(NumTrial), 
            Max = max(NumTrial)) %>% 
  ungroup %>%
  kable(digits=2)
## `summarise()` has grouped output by 'SN', 'cCong'. You can override using the `.groups` argument.
## `summarise()` has grouped output by 'cCong'. You can override using the `.groups` argument.
cCong cnRep Mean Median Min Max
v 1 55.06 56.0 48 58
v 2 56.00 56.0 51 59
v 3 56.44 58.0 50 59
v 4 55.97 56.0 44 60
v 0 56.19 56.5 49 60
nv 1 55.38 56.0 42 59
nv 2 56.22 57.0 51 60
nv 3 57.19 58.0 48 60
nv 4 57.22 57.0 52 60
nv 0 56.59 57.0 49 60

# check Distribution

# before trimming
den1 <- ggplot(ct1, aes(x=RT)) + 
  geom_density() + 
  theme_bw(base_size = 18) +
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank()) 

# after trimming
den2 <- ggplot(tt1, aes(x=RT)) + 
  geom_density() + 
  theme_bw(base_size = 18) + 
  labs(x = "Trimmed RT") + 
  theme(panel.grid.major = element_blank(), 
        panel.grid.minor = element_blank(),
        axis.text.y = element_blank(),
        axis.ticks.y = element_blank()) 
den1 + den2


ggplot(data=tt1, aes(x = SN, y = RT)) + 
  ggtitle("mean RT") +
  xlab("Participants") + ylab("RT") + theme_bw() +
  theme(text=element_text(size=14)) +
  theme(legend.title = element_text(face = 1,size = 15)) +
  geom_boxplot()



# subject-level, long format
t1rtL <- tt1 %>% filter(cnPair == 2) %>% group_by(SN, cCong, cnRep) %>%
  summarise(RT = mean(RT)) %>%
  ungroup()
## `summarise()` has grouped output by 'SN', 'cCong'. You can override using the `.groups` argument.
#t1rtL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t1rtW <- t1rtL %>% spread(key=cnRep, value = RT)
#t1rtW %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t1rtW1 <- t1rtL %>% # filter(cnRep == 4) %>% 
  spread(key=cCong, value = RT)
#t1rtW1 %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t1rtG <- t1rtL %>% group_by(cCong, cnRep) %>%
  summarise(RT.m = mean(RT), RT.sd = sd(RT)) %>%
  ungroup()
## `summarise()` has grouped output by 'cCong'. You can override using the `.groups` argument.
t1rtG$RT.se <- Rmisc::summarySEwithin(data = t1rtL, measurevar = "RT", 
                                      idvar = "SN", withinvars = c("cCong", "cnRep"))$se
t1rtG$RT.ci <- Rmisc::summarySEwithin(data = t1rtL, measurevar = "RT", 
                                      idvar = "SN", withinvars = c("cCong", "cnRep"))$ci
t1rtG <- t1rtG %>% 
  mutate(lower.ci = RT.m-RT.ci,
         upper.ci = RT.m+RT.ci)
t1rtG %>% dplyr::select(cCong, cnRep, RT.m) %>% 
  spread(key=cnRep, value=RT.m) %>% kable(digits=2)
cCong 1 2 3 4 0
v 758.18 704.12 667.56 678.37 733.66
nv 760.14 700.95 677.13 663.65 734.82


p1.L.rt.g <- t1rtG %>% filter(cnRep != 0)
p1.L.rt.g$Cong <- factor(p1.L.rt.g$cCong, labels=c("Violation", "Non-violation"))
p1.L.rt.g.plot1 <- ggplot(data=p1.L.rt.g , 
                              aes(x=cnRep, y=RT.m, ymin=RT.m-RT.ci, ymax=RT.m+RT.ci, color=Cong, shape=Cong)) +
  geom_point(size = 5, position = position_dodge(.3)) +
  geom_errorbar(width = .2, position = position_dodge(.3)) +
  geom_line(aes(group = Cong), position = position_dodge(.3), size=1) + 
  scale_color_manual(values = c("#ED7D31", "#70AD47")) + # c("red", "black")) +
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(500,1000), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  scale_x_discrete(labels = c("1", "2", "3", "4")) +
  labs(x = "Repetition", y = "Response Time (ms)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())

p1.L.rt.g <- t1rtG %>% filter(cnRep != 2, cnRep != 3, cnRep != 4)
p1.L.rt.g$Cong <- factor(p1.L.rt.g$cCong, labels=c("Violation", "Non-violation"))
p1.L.rt.g.plot2 <- ggplot(data=p1.L.rt.g , 
                          aes(x=cnRep, y=RT.m, ymin=RT.m-RT.ci, ymax=RT.m+RT.ci, color=Cong, shape=Cong)) +
  geom_point(size = 5, position = position_dodge(.3)) +
  geom_errorbar(width = .2, position = position_dodge(.3)) +
  geom_line(aes(group = Cong), position = position_dodge(.3), size=1) + 
  scale_color_manual(values = c("#ED7D31", "#70AD47")) + # c("red", "black")) +
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(500,1000), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  scale_x_discrete(labels = c("1st B", "B-Prime")) +
  labs(x = "Type", y = "Response Time (ms)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())
# p1.L.rt.g.plot1
# p1.L.rt.g.plot2

p1.rt.plot <- ggarrange(p1.L.rt.g.plot1, p1.L.rt.g.plot2, ncol = 2,
          labels=c("A) Categorization of B", "B) Categorization fo B/B`"))
p1.rt.plot


RM ANOVA 1 - 반복에 따른 priming 효과


t1rtL.1 <- t1rtL %>% filter(cnRep != 0)
t1.rt.aov1 <- aov_ez(id="SN", dv = "RT", data = t1rtL.1, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
# summary(t1.rt.aov1)
# anova(t1.rt.aov1, es = "pes")
nice(t1.rt.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 869.19 0.19 .006 .670
cnRep 2.26, 70.16 1856.46 77.77 *** .715 <.001
cCong:cnRep 2.32, 72.02 1261.63 1.70 .052 .185


RM ANOVA 1 - Post-hoc Test


### RM ANOVA 1 - post-hoc
t1.rt.m1 <- emmeans(t1.rt.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.rt.m1$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X1 -1.966 6.290 31 -0.313 0.757
v - nv X2 3.165 7.914 31 0.400 0.692
v - nv X3 -9.568 6.312 31 -1.516 0.140
v - nv X4 14.717 9.772 31 1.506 0.142
plot(t1.rt.m1, horizontal = T, comparisons = T)

t1.rt.m2 <- emmeans(t1.rt.aov1, pairwise ~ cnRep | cCong) # adjust="bon" , "tukey"
t1.rt.m2$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
X1 - X2 v 54.056 8.983 31 6.017 0.000
X1 - X3 v 90.618 9.704 31 9.338 0.000
X1 - X4 v 79.804 9.227 31 8.649 0.000
X2 - X3 v 36.563 8.987 31 4.068 0.002
X2 - X4 v 25.749 8.285 31 3.108 0.020
X3 - X4 v -10.814 7.597 31 -1.423 0.495
X1 - X2 nv 59.187 7.354 31 8.048 0.000
X1 - X3 nv 83.016 7.993 31 10.386 0.000
X1 - X4 nv 96.487 10.198 31 9.462 0.000
X2 - X3 nv 23.829 7.141 31 3.337 0.011
X2 - X4 nv 37.300 9.800 31 3.806 0.003
X3 - X4 nv 13.471 7.448 31 1.809 0.289
plot(t1.rt.m2, horizontal = T, comparisons = T)


RM ANOVA 2 - violation 여부에 따른 priming 효과 감소 (differentiation effect)


t1rtL.2 <- t1rtL %>% filter(cnRep !=1, cnRep != 2, cnRep !=0)
t1.rt.aov1 <- aov_ez(id="SN", dv = "RT", data = t1rtL.2, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
#summary(t1.rt.aov1)
#anova(t1.rt.aov1, es = "pes")
nice(t1.rt.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 969.52 0.22 .007 .643
cnRep 1, 31 615.18 0.09 .003 .764
cCong:cnRep 1, 31 1195.72 3.95 + .113 .056


RM ANOVA 2 - Post-hoc Test


### RM ANOVA 2 - post-hoc
t1.rt.m1 <- emmeans(t1.rt.aov1, pairwise ~ cnRep | cCong) # adjust="bon" , "tukey"
t1.rt.m1$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
X3 - X4 v -10.814 7.597 31 -1.423 0.165
X3 - X4 nv 13.471 7.448 31 1.809 0.080
plot(t1.rt.m1, horizontal = FALSE, comparisons = T)

t1.rt.m1 <- emmeans(t1.rt.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.rt.m1$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X3 -9.568 6.312 31 -1.516 0.140
v - nv X4 14.717 9.772 31 1.506 0.142
plot(t1.rt.m1, horizontal = FALSE, comparisons = T)


RM ANOVA 3 - general category priming 효과


t1rtL.3 <- t1rtL %>% filter(cnRep !=2, cnRep !=3, cnRep !=4)
t1.rt.aov1 <- aov_ez(id="SN", dv = "RT", data = t1rtL.3, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
# summary(t1.rt.aov1)
# anova(t1.rt.aov1, es = "pes")
nice(t1.rt.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 883.47 0.09 .003 .768
cnRep 1, 31 1031.18 19.27 *** .383 <.001
cCong:cnRep 1, 31 677.30 0.01 <.001 .931


RM ANOVA 3 - Post-hoc Test


t1.rt.m1 <- emmeans(t1.rt.aov1, pairwise ~ cnRep | cCong) # adjust="bon" , "tukey"
t1.rt.m1$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
X1 - X0 v 24.518 8.743 31 2.804 0.009
X1 - X0 nv 25.320 5.508 31 4.596 0.000
plot(t1.rt.m1, horizontal = FALSE, comparisons = T)

t1.rt.m2 <- emmeans(t1.rt.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.rt.m2$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X1 -1.966 6.290 31 -0.313 0.757
v - nv X0 -1.165 7.615 31 -0.153 0.879
plot(t1.rt.m2, horizontal = FALSE, comparisons = T)


2.1.2 RT analysis: Priming


t3.0 <- t1 %>% filter(Corr == 1) 

t3.1 <- t3.0 %>% filter(RT > 200 & RT < 10000) %>%
  group_by(SN) %>% # grouping by participants
  nest() %>%
  mutate(lbound = map(data, ~mean(.$RT)-3*sd(.$RT)),
         ubound = map(data, ~mean(.$RT)+3*sd(.$RT))) %>% # make new data (3sd cut)
  unnest(c(lbound, ubound))%>%
  unnest(data) %>%
  mutate(Outlier = (RT < lbound)|(RT > ubound)) %>% # set outlier
  filter(Outlier == FALSE) %>% # filtering outlier
  ungroup()

# 100-100*(nrow(t3.1)/nrow(t3.0))

t3.1 <- t3.1 %>% filter(cnRep!=0, cnRep!=1, cnRep!=2)
glimpse(t3.1, width = 50)
## Rows: 10,873
## Columns: 18
## $ SN      <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ Trial   <int> 25, 26, 27, 29, 30, 31, 32, 34, …
## $ Block   <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ bTrial  <int> 25, 26, 27, 29, 30, 31, 32, 34, …
## $ cCong   <fct> nv, nv, v, v, v, nv, nv, nv, nv,…
## $ cnRep   <fct> 3, 3, 3, 3, 3, 3, 3, 3, 3, 3, 4,…
## $ cnPair  <fct> 1, 2, 1, 1, 2, 1, 2, 1, 2, 2, 2,…
## $ cIdt    <fct> 4, 4, 2, 5, 5, 3, 3, 6, 6, 1, 4,…
## $ cCat    <fct> 2, 1, 1, 2, 1, 2, 1, 1, 2, 2, 1,…
## $ Im_nExm <int> 3, 1, 3, 3, 1, 3, 1, 3, 1, 1, 1,…
## $ Im_Idx  <fct> 209, 111, 47, 181, 72, 157, 127,…
## $ Im_Name <fct> 209_3_out.jpg, 111_1_in.jpg, 047…
## $ Resp    <int> 2, 1, 1, 2, 1, 2, 1, 1, 2, 2, 1,…
## $ RT      <dbl> 589.8, 546.3, 935.2, 591.1, 525.…
## $ Corr    <dbl> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ lbound  <dbl> 199.812, 199.812, 199.812, 199.8…
## $ ubound  <dbl> 996.88, 996.88, 996.88, 996.88, …
## $ Outlier <lgl> FALSE, FALSE, FALSE, FALSE, FALS…

t3.2 <- t3.1 %>% filter(cnPair == 2) %>% 
  dplyr::select(SN, cCong, cnRep, cIdt, Im_Name, RT)

t3.3 <- t3.2 %>% group_by(SN) %>% 
  spread(key = cnRep, value = RT) %>% ungroup()

t3.4 <- t3.3 %>% group_by(SN, cCong, cIdt) %>% mutate(prime = `4` - `3`) %>% ungroup()
t3.4 <- na.omit(t3.4)
glimpse(t3.4)
## Rows: 3,505
## Columns: 7
## $ SN      <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1,…
## $ cCong   <fct> v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v,…
## $ cIdt    <fct> 1, 5, 8, 11, 12, 13, 14, 15, 16, 21, 22, 26, 27, 30, 31, 33, 3…
## $ Im_Name <fct> 193_1_out.jpg, 072_1_in.jpg, 124_1_in.jpg, 150_1_out.jpg, 212_…
## $ `3`     <dbl> 520.2, 525.2, 755.2, 722.3, 487.2, 515.7, 552.4, 518.7, 535.0,…
## $ `4`     <dbl> 471.9, 560.4, 588.6, 549.8, 568.8, 513.3, 536.6, 551.9, 522.2,…
## $ prime   <dbl> -48.3, 35.2, -166.6, -172.5, 81.6, -2.4, -15.8, 33.2, -12.8, -…

# subject-level, long format
t3rtL <- t3.4 %>% group_by(SN, cCong) %>%
  summarise(prime = mean(prime)) %>%
  ungroup()
## `summarise()` has grouped output by 'SN'. You can override using the `.groups` argument.
#t3rtL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t3rtW <- t3rtL %>% spread(key=cCong, value = prime)
#t3rtW %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t3rtG <- t3rtL %>% group_by(cCong) %>%
  summarise(prime.m = mean(prime), prime.sd = sd(prime)) %>%
  ungroup()
t3rtG$prime.se <- Rmisc::summarySEwithin(data = t3rtL, measurevar = "prime", 
                                         idvar = "SN", withinvars = c("cCong"))$se
t3rtG$prime.ci <- Rmisc::summarySEwithin(data = t3rtL, measurevar = "prime", 
                                         idvar = "SN", withinvars = c("cCong"))$ci
t3rtG <- t3rtG %>% 
  mutate(lower.ci = prime.m-prime.ci,
         upper.ci = prime.m+prime.ci)
t3rtG %>% kable(digits=2)
cCong prime.m prime.sd prime.se prime.ci lower.ci upper.ci
v 10.68 42.24 8.41 17.16 -6.48 27.84
nv -16.03 41.28 8.41 17.16 -33.19 1.13



# plot
p1.L.rt3.g <- t3rtG
p1.L.rt3.g$Cong <- factor(p1.L.rt3.g$cCong, labels=c("Violation", "Non-violation"))
p1.L.rt3.l <- t3rtL
p1.L.rt3.l$Cong <- factor(p1.L.rt3.l$cCong, labels=c("Violation", "Non-violation"))
p1.L.rt3.w <- p1.L.rt3.l %>% dplyr::select(SN, Cong, prime) %>% spread(key =Cong, value=prime)

p1.L.rt3.g.plot3 <- ggplot() +
  geom_hline(yintercept=0, linetype='dashed', color='gray60', alpha =1, size=1) +
  geom_point(data=p1.L.rt3.l, aes(x=Cong, y=prime), size = 1, position = position_dodge(.8),  color="gray80") +
  geom_segment(data=p1.L.rt3.w, inherit.aes = FALSE,
               aes(x=1, y=`Violation`,
                   xend=2, yend=`Non-violation`),
               color="gray80") +
  geom_pointrange(data=p1.L.rt3.g, aes(x = Cong, y=prime.m, ymin = prime.m-prime.ci, ymax = prime.m+prime.ci, color = Cong),
                  position = position_dodge(0.80), size = 1.5, show.legend = FALSE) +
  scale_color_manual(values = c("#ED7D31", "#70AD47")) + 
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(-150,+150), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  scale_x_discrete(labels = c("Violation", "Non-violation")) +
  labs(x = "Type", y = "4th - 3rd (ms)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())

p1.L.rt3.g.plot3

# p1.rt.plot2 <- ggarrange(p1.L.rt.g.plot1, p1.L.rt.g.plot2,p1.L.rt3.g.plot3, ncol = 3,
#                          labels=c("A) Categorization of B", "B) Categorization fo B/B`", "C) Differentiation Effect"))
# p1.rt.plot2


p1.L.rt3.g.plot3.sn <- ggplot() +
  geom_hline(yintercept=0,  color='gray60', alpha =1, size=1) +
  stat_summary(data=p1.L.rt3.l, aes(x=SN, y=prime, fill=Cong), 
               fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02) + # , show.legend = FALSE  colour="black", 
  
  scale_fill_manual(values = c("#ED7D31", "#70AD47"),
                    labels = c("Violation", "Non-violation")) + 
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(-120,+120), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  # scale_x_discrete(labels = c("Violation", "Non-violation")) +
  labs(x = "Type", y = "4th - 3rd (ms)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = 'top', # c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())

p1.L.rt3.g.plot3.sn


Paired T-Test 1 - violation 여부에 따른 priming 약화 효과 (differentiation effect)


t3rtL.1 <- t3rtL
t3.rt.ttest <- t.test(prime ~ cCong, data=t3rtL.1,
                      alternative = c("two.sided"),
                      paired=T, var.equal=F,
                      conf.level=.95)
t3.rt.ttest
## 
##  Paired t-test
## 
## data:  prime by cCong
## t = 2.2441, df = 31, p-value = 0.03211
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##   2.43423 50.97299
## sample estimates:
## mean of the differences 
##                26.70361
t3.rt.rmanova <- aov_ez(id="SN", dv = "prime", data = t3rtL.1, within = c("cCong"))
nice(t3.rt.rmanova, es ="pes")
## Anova Table (Type 3 tests)
## 
## Response: prime
##   Effect    df     MSE      F  pes p.value
## 1  cCong 1, 31 2265.60 5.04 * .140    .032
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '+' 0.1 ' ' 1
t3.rt.m1 <- emmeans(t3.rt.rmanova, pairwise ~ cCong, adjust="tukey") # adjust="bon" , "tukey"
t3.rt.m1$contrasts %>% kable(digits=3)
contrast estimate SE df t.ratio p.value
v - nv 26.704 11.9 31 2.244 0.032
plot(t3.rt.m1, horizontal = FALSE, comparisons = T)


2.1.3 RT analysis: Priming with B’ RT



t4.0 <- t1 %>% filter(Corr == 1)

t4.1 <- t4.0 %>% filter(RT > 200 & RT < 10000) %>%
  group_by(SN) %>% # grouping by participants
  nest() %>%
  mutate(lbound = map(data, ~mean(.$RT)-3*sd(.$RT)),
         ubound = map(data, ~mean(.$RT)+3*sd(.$RT))) %>% # make new data (3sd cut)
  unnest(c(lbound, ubound))%>%
  unnest(data) %>%
  mutate(Outlier = (RT < lbound)|(RT > ubound)) %>% # set outlier
  filter(Outlier == FALSE) %>% # filtering outlier
  ungroup()

100-100*(nrow(t4.1)/nrow(t4.0))
## [1] 1.812514

t4.1 <- t4.1 %>% filter(cnRep!=1, cnRep!=2)
# t4.1$cnRep <- as.double(t4.1$cnRep)
t4.2 <- t4.1 %>% filter(cnPair == 2) %>% 
  dplyr::select(SN, cCong, cnRep, cIdt, Im_Idx, RT)

t4.3 <- t4.2 %>% group_by(SN) %>% 
  spread(key = cnRep, value = RT) %>% ungroup()


t4.4 <- t4.3 %>% group_by(SN, cIdt) %>% mutate(prime = `4` - `3`) %>% ungroup()
t4.4 <- na.omit(t4.4)
t4.4 <- t4.4 %>% dplyr::select (SN, cCong, cIdt, Im_Idx, `0`, prime)
glimpse(t4.4)
## Rows: 3,345
## Columns: 6
## $ SN     <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
## $ cCong  <fct> v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, v, …
## $ cIdt   <fct> 1, 5, 11, 12, 13, 14, 16, 21, 22, 26, 27, 30, 31, 33, 34, 35, 3…
## $ Im_Idx <fct> 193, 72, 150, 212, 234, 220, 122, 149, 109, 180, 207, 95, 80, 7…
## $ `0`    <dbl> 663.1, 631.0, 596.8, 548.1, 534.3, 752.2, 516.5, 685.6, 572.1, …
## $ prime  <dbl> -48.3, 35.2, -172.5, 81.6, -2.4, -15.8, -12.8, -92.9, -21.5, -6…

# subject-level, long format
t4rtL <- t4.4 %>% group_by(SN, cCong) %>%
  summarise(prime = mean(prime), bpr = mean(`0`)) %>%
  ungroup()
## `summarise()` has grouped output by 'SN'. You can override using the `.groups` argument.
#t4rtL %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t4rtG <- t4rtL %>% group_by(cCong) %>%
  summarise(prime.m = mean(prime), prime.sd = sd(prime), bpr.m = mean(bpr), brp.sd =sd(bpr)) %>%
  ungroup()
t3rtG %>% kable(digits=2)
cCong prime.m prime.sd prime.se prime.ci lower.ci upper.ci
v 10.68 42.24 8.41 17.16 -6.48 27.84
nv -16.03 41.28 8.41 17.16 -33.19 1.13



# correlation by Subject
t4rtL.1 <- t4rtL %>% filter(cCong == 'v')
cor.test(formula = ~ prime + bpr, data = t4rtL.1,
         method = "pearson", alternative = "two.sided")
## 
##  Pearson's product-moment correlation
## 
## data:  prime and bpr
## t = 1.0021, df = 30, p-value = 0.3243
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1800267  0.4974426
## sample estimates:
##       cor 
## 0.1799644



t4rtL.2 <- t4rtL %>% filter(cCong == 'nv')
cor.test(formula = ~ prime + bpr, data = t4rtL.2,
         method = "pearson", alternative = "two.sided")
## 
##  Pearson's product-moment correlation
## 
## data:  prime and bpr
## t = 0.7508, df = 30, p-value = 0.4586
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.2234696  0.4625946
## sample estimates:
##       cor 
## 0.1358068



p1.rt.corr.plot1 <- ggarrange(p1.L.rt.corr.plot1, p1.L.rt.corr.plot2, ncol = 2,
                         labels=c("A) Violation Condition", "B) Non-violation Condition"))
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
p1.rt.corr.plot1


2.1.4 RT analysis: correlation, aggregation

t4.4.1 <- t4.4 %>% filter(cCong =="v")
cor.test(formula = ~ prime + `0`, data = t4.4.1,
         method = "pearson", alternative = "two.sided")
## 
##  Pearson's product-moment correlation
## 
## data:  prime and 0
## t = 0.63219, df = 1646, p-value = 0.5273
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.03273072  0.06381923
## sample estimates:
##        cor 
## 0.01558057
ggplot(t4.4.1, aes(x=`0`, y=prime)) +
  geom_point(size = 4) +
  geom_smooth(method=lm) +
  labs(x = "b` RT", 
       y = "4th - 3rd B RT") +
  theme_bw(base_size = 18) +
  theme(panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        aspect.ratio = 1)
## `geom_smooth()` using formula 'y ~ x'

t4.4.2 <- t4.4 %>% filter(cCong =="nv")
cor.test(formula = ~ prime + `0`, data = t4.4.2,
         method = "pearson", alternative = "two.sided")
## 
##  Pearson's product-moment correlation
## 
## data:  prime and 0
## t = 2.0711, df = 1695, p-value = 0.0385
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.002663411 0.097592241
## sample estimates:
##       cor 
## 0.0502413
ggplot(t4.4.2, aes(x=`0`, y=prime)) +
  geom_point(size = 4) +
  geom_smooth(method=lm) +
  labs(x = "b` RT", 
       y = "4th - 3rd B RT") +
  theme_bw(base_size = 18) +
  theme(panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        aspect.ratio = 1)
## `geom_smooth()` using formula 'y ~ x'


2.1.5 Accuracy analysis


t1acc.sn <- t1 %>% filter(cnPair == 2) %>% group_by(SN, cCong) %>%
  dplyr::summarise(Acc = mean(Corr)*100) %>%
  ungroup() # %>% spread(key=cCong,value=Acc)
## `summarise()` has grouped output by 'SN'. You can override using the `.groups` argument.
#


ggplot(data=t1acc.sn, aes(x = SN, y = Acc)) + 
  ggtitle("mean Accuracy") +
  xlab("Participant") + ylab("Acc") + theme_bw() +
  theme(text=element_text(size=14)) +
  theme(legend.title = element_text(face = 1,size = 15)) +
  geom_boxplot()


t1accL <- t1 %>% filter(cnPair == 2) %>% 
  group_by(SN, cCong, cnRep) %>%
  dplyr::summarise(Acc = mean(Corr)*100) %>%
  ungroup()
## `summarise()` has grouped output by 'SN', 'cCong'. You can override using the `.groups` argument.
# t1accL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t1accW1 <- t1accL %>% spread(key=cCong, value = Acc)
# t1accW1 %>% kable(digits=2)

t1accW2 <- t1accL %>% spread(key=cnRep, value = Acc)
# t1accW2 %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t1accG <- t1accL %>% group_by(cCong, cnRep) %>%
  summarise(Acc.m = mean(Acc), Acc.sd = sd(Acc)) %>%
  ungroup()
## `summarise()` has grouped output by 'cCong'. You can override using the `.groups` argument.
t1accG$Acc.se <- Rmisc::summarySEwithin(data = t1accL, measurevar = "Acc", 
                                        idvar = "SN", withinvars = c("cCong", "cnRep"))$se
t1accG$Acc.ci <- Rmisc::summarySEwithin(data = t1accL, measurevar = "Acc", 
                                        idvar = "SN", withinvars = c("cCong", "cnRep"))$ci
t1accG <- t1accG %>% 
  mutate(lower.ci = Acc.m-Acc.ci,
         upper.ci = Acc.m+Acc.ci)
t1accG %>% 
  dplyr::select(cCong, cnRep, Acc.m) %>% 
  spread(key=cnRep, value= Acc.m) %>% kable(digits=2)
cCong 1 2 3 4 0
v 93.85 95.16 95.42 95.21 95.57
nv 94.84 95.42 96.35 96.51 96.30

<br

p1.L.acc.g <- t1accG %>% filter(cnRep != 0)
p1.L.acc.g$Cong <- factor(p1.L.acc.g$cCong, labels=c("Violation", "Non-violation"))
p1.L.acc.g.plot1 <- ggplot(data=p1.L.acc.g , 
                           aes(x=cnRep, y=Acc.m, ymin=Acc.m-Acc.ci, ymax=Acc.m+Acc.ci, color=Cong, shape=Cong)) +
  geom_point(size = 5, position = position_dodge(.3)) +
  geom_errorbar(width = .2, position = position_dodge(.3)) +
  geom_line(aes(group = Cong), position = position_dodge(.3), size=1) + 
  scale_color_manual(values = c("#ED7D31", "#70AD47")) + # c("red", "black")) +
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(80,110), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  scale_x_discrete(labels = c("1", "2", "3", "4")) +
  labs(x = "Repetition", y = "Accuracy (%)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())

p1.L.acc.g <- t1accG %>% filter(cnRep != 2, cnRep != 3, cnRep != 4)
p1.L.acc.g$Cong <- factor(p1.L.acc.g$cCong, labels=c("Violation", "Non-violation"))
p1.L.acc.g.plot2 <- ggplot(data=p1.L.acc.g , 
                           aes(x=cnRep, y=Acc.m, ymin=Acc.m-Acc.ci, ymax=Acc.m+Acc.ci, color=Cong, shape=Cong)) +
  geom_point(size = 5, position = position_dodge(.3)) +
  geom_errorbar(width = .2, position = position_dodge(.3)) +
  geom_line(aes(group = Cong), position = position_dodge(.3), size=1) + 
  scale_color_manual(values = c("#ED7D31", "#70AD47")) + # c("red", "black")) +
  # scale_shape_manual(values = c("red", "black"), labels = c("Experimental", "Control")) +
  coord_cartesian(ylim = c(80,110), clip = "on") +
  # scale_y_continuous(breaks = seq(0,80)) +
  scale_x_discrete(labels = c("1st B", "B-Prime")) +
  labs(x = "Type", y = "Accuracy (%)") +
  theme_bw(base_size = 18) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        legend.position = c(0.65, 0.85),
        legend.title = element_blank(),
        # aspect.ratio = 0.6,
        legend.background = element_blank(),
        plot.margin = margin(1, 0.3, 0.3, 0.5, "cm"),
        legend.key = element_blank())

p1.acc.plot <- ggarrange(p1.L.acc.g.plot1, p1.L.acc.g.plot2, ncol = 2,
                         labels=c("A) Categorization of B", "B) Categorization fo B/B`"))
p1.acc.plot


RM ANOVA


t1accL.1 <- t1accL # %>% filter(SN!=4) 
t1.acc.aov1 <- aov_ez(id="SN", dv = "Acc", data = t1accL.1, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
#summary(t1.acc.aov1)
#anova(t1.acc.aov1, es = "pes")
nice(t1.acc.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 14.98 3.80 + .109 .060
cnRep 3.09, 95.74 7.73 4.91 ** .137 .003
cCong:cnRep 3.16, 97.87 6.57 0.46 .015 .722


RM ANOVA 1 - post-hoc


t1.acc.m1 <- emmeans(t1.acc.aov1, pairwise ~ cnRep | cCong) # adjust="bon" , "tukey"
t1.acc.m1$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
X1 - X2 v -1.302 0.634 31 -2.054 0.265
X1 - X3 v -1.562 0.617 31 -2.533 0.109
X1 - X4 v -1.354 0.582 31 -2.328 0.163
X1 - X0 v -1.719 0.615 31 -2.796 0.062
X2 - X3 v -0.260 0.627 31 -0.416 0.993
X2 - X4 v -0.052 0.663 31 -0.079 1.000
X2 - X0 v -0.417 0.594 31 -0.701 0.955
X3 - X4 v 0.208 0.593 31 0.351 0.997
X3 - X0 v -0.156 0.542 31 -0.289 0.998
X4 - X0 v -0.365 0.528 31 -0.691 0.957
X1 - X2 nv -0.573 0.659 31 -0.869 0.906
X1 - X3 nv -1.510 0.542 31 -2.789 0.063
X1 - X4 nv -1.667 0.698 31 -2.388 0.146
X1 - X0 nv -1.458 0.701 31 -2.080 0.254
X2 - X3 nv -0.937 0.449 31 -2.090 0.250
X2 - X4 nv -1.094 0.583 31 -1.877 0.350
X2 - X0 nv -0.885 0.467 31 -1.895 0.341
X3 - X4 nv -0.156 0.493 31 -0.317 0.998
X3 - X0 nv 0.052 0.505 31 0.103 1.000
X4 - X0 nv 0.208 0.634 31 0.329 0.997
plot(t1.acc.m1, horizontal = FALSE, comparisons = T)

t1.acc.m2 <- emmeans(t1.acc.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.acc.m2$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X1 -0.990 0.681 31 -1.453 0.156
v - nv X2 -0.260 0.508 31 -0.512 0.612
v - nv X3 -0.937 0.639 31 -1.467 0.152
v - nv X4 -1.302 0.829 31 -1.571 0.126
v - nv X0 -0.729 0.643 31 -1.133 0.266
plot(t1.acc.m2, horizontal = FALSE, comparisons = T)


RM-ANOVA 2


t1accL.2 <- t1accL %>% filter(cnRep !=1, cnRep != 2, cnRep !=0)
t1.acc.aov1 <- aov_ez(id="SN", dv = "Acc", data = t1accL.2, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
# summary(t1.acc.aov1)
# anova(t1.acc.aov1, es = "pes")
nice(t1.acc.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 11.99 3.35 + .097 .077
cnRep 1, 31 3.96 0.01 <.001 .941
cCong:cnRep 1, 31 5.54 0.19 .006 .664


RM ANOVA 2 - post-hoc


t1.acc.m1 <- emmeans(t1.acc.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.acc.m1$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X3 -0.937 0.639 31 -1.467 0.152
v - nv X4 -1.302 0.829 31 -1.571 0.126
plot(t1.acc.m1, horizontal = FALSE, comparisons = T)

t1.acc.m1 <- emmeans(t1.acc.aov1, pairwise ~ cnRep) # adjust="bon" , "tukey"
t1.acc.m1$contrasts %>% kable(digits=3)
contrast estimate SE df t.ratio p.value
X3 - X4 0.026 0.352 31 0.074 0.941
plot(t1.acc.m1, horizontal = FALSE, comparisons = T)


Paired T-Test


t1accL.3 <- t1accL %>% filter(cnRep == 0)
t1.acc.ttest <- t.test(Acc ~ cCong, data=t1accL.3,
                       alternative = c("two.sided"),
                       paired=T, var.equal=F,
                       conf.level=.95)
t1.acc.ttest 
## 
##  Paired t-test
## 
## data:  Acc by cCong
## t = -1.1331, df = 31, p-value = 0.2658
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -2.0415686  0.5832353
## sample estimates:
## mean of the differences 
##              -0.7291667


RM-ANOVA 3


t1accL.3 <- t1accL %>% filter(cnRep !=2, cnRep !=3, cnRep !=4)
t1.acc.aov1 <- aov_ez(id="SN", dv = "Acc", data = t1accL.3, within = c("cCong", "cnRep"))
# test_sphericity(t2.acc.aov1)
# summary(t1.acc.aov1)
# anova(t1.acc.aov1, es = "pes")
nice(t1.acc.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 9.79 2.41 .072 .130
cnRep 1, 31 9.65 8.37 ** .213 .007
cCong:cnRep 1, 31 4.26 0.13 .004 .724


RM-ANOVA 3 - post hoc


### RM ANOVA 3 - post-hoc
t1.acc.m1 <- emmeans(t1.acc.aov1, pairwise ~ cnRep | cCong) # adjust="bon" , "tukey"
t1.acc.m1$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
X1 - X0 v -1.719 0.615 31 -2.796 0.009
X1 - X0 nv -1.458 0.701 31 -2.080 0.046
plot(t1.acc.m1, horizontal = FALSE, comparisons = T)

t1.rt.m2 <- emmeans(t1.rt.aov1, pairwise ~ cCong | cnRep) # adjust="bon" , "tukey"
t1.rt.m2$contrasts %>% kable(digits=3)
contrast cnRep estimate SE df t.ratio p.value
v - nv X1 -1.966 6.290 31 -0.313 0.757
v - nv X0 -1.165 7.615 31 -0.153 0.879
plot(t1.rt.m2, horizontal = FALSE, comparisons = T)




2.2 Phase 2 - Recognition Test


t2 <- read.csv("data/reactMEM_p2_all.csv", header = T)

glimpse(t2, width = 70)
## Rows: 7,680
## Columns: 13
## $ SN       <int> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
## $ Trial    <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, …
## $ cType    <int> 2, 1, 1, 2, 1, 2, 2, 1, 2, 2, 1, 2, 2, 1, 2, 1, 2, …
## $ cCong    <int> 999, 1, 2, 999, 1, 999, 999, 2, 999, 999, 2, 999, 9…
## $ cIdt     <int> 999, 117, 66, 999, 60, 999, 999, 95, 999, 999, 113,…
## $ cCat     <int> 2, 2, 1, 1, 2, 2, 2, 1, 2, 1, 2, 2, 1, 1, 2, 1, 2, …
## $ xIm_nExm <int> 999, 2, 2, 999, 2, 999, 999, 2, 999, 999, 2, 999, 9…
## $ xIm_Idx  <int> 241, 227, 71, 301, 191, 341, 303, 27, 252, 282, 143…
## $ dIm_Name <chr> "241_1_out.jpg", "227_2_out.jpg", "071_2_in.jpg", "…
## $ Resp     <int> 1, 4, 2, 2, 2, 1, 1, 4, 1, 2, 4, 2, 1, 3, 2, 2, 1, …
## $ RT       <dbl> 2.8079, 2.2212, 2.0978, 1.2932, 1.6856, 1.5962, 1.3…
## $ Corr     <int> 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, …
## $ Conf     <int> 2, 2, 1, 1, 1, 2, 2, 2, 2, 1, 2, 1, 2, 1, 1, 1, 2, …

# ------- Variables ------- #
# $ SN       <int> 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99, 99…
# $ Trial    <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15,…
# $ cType    <int> 1, 1, 1, 2, 2, 1, 2, 2, 1, 1, 1, 1, 2, 2, 2, 2, 1,…
# $ cCong    <int> 1, 2, 1, 999, 999, 1, 999, 999, 2, 1, 1, 1, 999, 9…
# $ cIdt     <int> 3, 105, 43, 999, 999, 89, 999, 999, 103, 81, 23, 2…
# $ cCat     <int> 2, 1, 2, 2, 2, 1, 2, 1, 2, 1, 1, 1, 1, 1, 2, 2, 1,…
# $ xIm_nExm <int> 2, 2, 2, 999, 999, 2, 999, 999, 2, 2, 2, 2, 999, 9…
# $ xIm_Idx  <int> 199, 130, 197, 271, 259, 4, 278, 274, 192, 73, 62,…
# $ dIm_Name <chr> " 199_2_out.jpg", " 130_2_in.jpg", " 197_2_out.jpg…
# $ Resp     <int> 2, 4, 3, 1, 1, 4, 4, 1, 2, 4, 2, 1, 1, 1, 1, 3, 4,…
# $ RT       <dbl> 62.2271, 2.4708, 2.2955, 5.4393, 1.4200, 1.1727, 1…
# $ Corr     <int> 0, 1, 1, 1, 1, 1, 0, 1, 0, 1, 0, 0, 1, 1, 1, 0, 1,…
# $ Conf     <int> 1, 2, 1, 2, 2, 2, 2, 2, 1, 2, 1, 2, 2, 2, 2, 1, 2,…

length(unique(t2$SN))
## [1] 32

# check number of trials for each condition/SN
table(t2$SN)
## 
##   1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19  20 
## 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 240 
##  21  22  23  24  25  26  27  28  29  30  31  32 
## 240 240 240 240 240 240 240 240 240 240 240 240
table(t2$cType, t2$SN) 
##    
##       1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18  19
##   1 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120
##   2 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120
##    
##      20  21  22  23  24  25  26  27  28  29  30  31  32
##   1 120 120 120 120 120 120 120 120 120 120 120 120 120
##   2 120 120 120 120 120 120 120 120 120 120 120 120 120
table(t2$cCong, t2$SN) 
##      
##         1   2   3   4   5   6   7   8   9  10  11  12  13  14  15  16  17  18
##   1    60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60
##   2    60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60  60
##   999 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120 120
##      
##        19  20  21  22  23  24  25  26  27  28  29  30  31  32
##   1    60  60  60  60  60  60  60  60  60  60  60  60  60  60
##   2    60  60  60  60  60  60  60  60  60  60  60  60  60  60
##   999 120 120 120 120 120 120 120 120 120 120 120 120 120 120
# table(t2$cCat, t2$SN) 

# change class of main factors: double to factor
t2$SN = factor(t2$SN)
t2$cType = factor(t2$cType, levels=c(1,2), labels=c("old","new"))
t2$cCong = factor(t2$cCong, levels=c(1,2,999), labels=c("v","nv","new")) 
t2$cIdt = factor(t2$cIdt)
# t2$cCat = factor(t2$cCat)
t2$dIm_Name = factor(t2$dIm_Name)
t2$Resp <- factor(t2$Resp);
t2$mResp <- factor(t2$Resp, levels=c(1,2,3,4), labels=c("new", "new", "old", "old"))
t2$RT <- t2$RT*1000;
t2$Corr <- as.numeric(t2$Corr==1)
t2$Conf <- as.double(t2$Conf)
# t2$Corr[t2$Conf==2] <- 1
# t2$Corr[t2$cType=="old" & t2$mResp=="old" &  t2$Conf==2] <- 1
# t2$Corr[t2$cType=="old" & t2$mResp=="old" & t2$Conf==1] <- 0

t2$p_resp = factor(t2$cType, levels=c('old','new'), labels=c("p.old","p.new"))
t2$r_resp = factor(t2$mResp, levels=c('old','new'), labels=c("r.old", "r.new"))

glimpse(t2)
## Rows: 7,680
## Columns: 16
## $ SN       <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1…
## $ Trial    <int> 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18…
## $ cType    <fct> new, old, old, new, old, new, new, old, new, new, old, new, n…
## $ cCong    <fct> new, v, nv, new, v, new, new, nv, new, new, nv, new, new, nv,…
## $ cIdt     <fct> 999, 117, 66, 999, 60, 999, 999, 95, 999, 999, 113, 999, 999,…
## $ cCat     <int> 2, 2, 1, 1, 2, 2, 2, 1, 2, 1, 2, 2, 1, 1, 2, 1, 2, 2, 2, 2, 2…
## $ xIm_nExm <int> 999, 2, 2, 999, 2, 999, 999, 2, 999, 999, 2, 999, 999, 2, 999…
## $ xIm_Idx  <int> 241, 227, 71, 301, 191, 341, 303, 27, 252, 282, 143, 280, 245…
## $ dIm_Name <fct> 241_1_out.jpg, 227_2_out.jpg, 071_2_in.jpg, 301_1_in.jpg, 191…
## $ Resp     <fct> 1, 4, 2, 2, 2, 1, 1, 4, 1, 2, 4, 2, 1, 3, 2, 2, 1, 2, 4, 2, 2…
## $ RT       <dbl> 2807.9, 2221.2, 2097.8, 1293.2, 1685.6, 1596.2, 1394.6, 1831.…
## $ Corr     <dbl> 1, 1, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 1, 1, 1, 1, 1…
## $ Conf     <dbl> 2, 2, 1, 1, 1, 2, 2, 2, 2, 1, 2, 1, 2, 1, 1, 1, 2, 1, 2, 1, 1…
## $ mResp    <fct> new, old, new, new, new, new, new, old, new, new, old, new, n…
## $ p_resp   <fct> p.new, p.old, p.old, p.new, p.old, p.new, p.new, p.old, p.new…
## $ r_resp   <fct> r.new, r.old, r.new, r.new, r.new, r.new, r.new, r.old, r.new…


2.2.1 Accuracy Analysis: Raw Data


t2acc.sn <- t2 %>% group_by(SN, cType, cCong) %>%
  dplyr::summarise(Acc = mean(Corr)*100) %>%
  ungroup() %>% spread(key=cType,value=Acc)
## `summarise()` has grouped output by 'SN', 'cType'. You can override using the `.groups` argument.
# t2acc.sn %>% kable(digits=2)

t2acc.snL <- t2 %>% 
  filter(cType=='old') %>%
  group_by(SN, cCong) %>%
  dplyr::summarise(Acc = mean(Corr)*100) %>%
  ungroup()
## `summarise()` has grouped output by 'SN'. You can override using the `.groups` argument.
#t2acc.snL %>% kable(digits=2)

ggplot(data=t2acc.snL, aes(x = SN, y = Acc)) + 
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02) + # , show.legend = FALSE  colour="black",
  geom_hline(yintercept=50, linetype='dashed', color='darkred', alpha =1, size=1) +
  coord_cartesian(ylim = c(0, 100), clip = "on") +
  ggtitle("mean Accuracy (studied)") +
  xlab("Participants") + ylab("Acc") + theme_bw() +
  theme(text=element_text(size=14)) +
  theme(legend.title = element_text(face = 1,size = 15))


# subject-level, long format (SN/Btw/Block)
t2accL <- t2 %>% group_by(SN, cType, cCong) %>%
  dplyr::summarise(Acc = mean(Corr)*100) %>%
  ungroup()
## `summarise()` has grouped output by 'SN', 'cType'. You can override using the `.groups` argument.
# t2accL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t2accW1 <- t2accL %>% spread(key=cCong, value = Acc)
# t2accW1 %>% kable(digits=2)

t2accW2 <- t2accL %>% spread(key=cType, value = Acc)
# t2accW2 %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t2accG <- t2accL %>% group_by(cType, cCong) %>%
  summarise(Acc.m = mean(Acc), Acc.sd = sd(Acc)) %>%
  ungroup()
## `summarise()` has grouped output by 'cType'. You can override using the `.groups` argument.
t2accG$Acc.se <- Rmisc::summarySEwithin(data = t2accL, measurevar = "Acc", 
                                        idvar = "SN", withinvars = c("cType", "cCong"))$se
t2accG$Acc.ci <- Rmisc::summarySEwithin(data = t2accL, measurevar = "Acc", 
                                        idvar = "SN", withinvars = c("cType", "cCong"))$ci
t2accG <- t2accG %>% 
  mutate(lower.ci = Acc.m-Acc.ci,
         upper.ci = Acc.m+Acc.ci)
# for between-subject design. check help(summarySE)
t2accG %>% kable(digits=2)
cType cCong Acc.m Acc.sd Acc.se Acc.ci lower.ci upper.ci
old v 56.15 15.71 1.41 2.88 53.26 59.03
old nv 52.40 16.34 1.51 3.09 49.31 55.48
new new 94.45 6.25 2.51 5.12 89.34 99.57


t2.plot1 <- ggplot(data=t2accL, aes(x=cType, y=Acc, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02) + # , show.legend = FALSE  colour="black", 
  geom_point(data=t2accL, aes(x=cType, y=Acc, fill=cCong), 
             position=position_dodge(width = 0.8), 
             size=2, show.legend = FALSE, color="gray60") + # , 
  geom_segment(data=filter(t2accW1, cType == "old"), inherit.aes = FALSE,
               aes(x=0.8, y=filter(t2accW1, cType == "old")$v,
                   xend=1.2, yend=filter(t2accW1, cType == "old")$nv),
               color="gray60") +
  geom_segment(data=filter(t2accW1, cType == "new"), inherit.aes = FALSE,
               aes(x=1.8, y=filter(t2accW1, cType == "new")$v,
                   xend=2.2, yend=filter(t2accW1, cType == "new")$nv),
               color="gray60") +
  geom_pointrange(data=t2accG, aes(x = cType, y=Acc.m, ymin = Acc.m-Acc.ci, ymax = Acc.m+Acc.ci),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) +
  scale_x_discrete(labels=c("Studied","Lure")) +
  scale_fill_manual(values = c("#D2E6BD", "#4AA7B4", "#235796"),# c("#6798D4", "#D87B49", "#4FA77C"), 
                    labels = c("Violation", "Non-violation", "Lure")) +
  # scale_fill_manual(values = c("#feb24c", "#91bfdb"), #6798D4 #D87B49 #4FA77C
                    # labels = c("Violation", "Non-violation")) +
  coord_cartesian(ylim = c(0, 120), clip = "on") +
  scale_y_continuous(breaks = c(0,25,50,75,100)) +
  labs(x = "Stimuli Type", y = "Accuracy (%)", fill ="Violation Condition") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position="top")
        # legend.position=c(0.8, 0.90))
t2.plot1


t2accL.sn <- t2accL %>% filter(cCong!="new")
t2.plot1.sn <- ggplot(data=t2accL.sn, aes(x=SN, y=Acc, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02) + # , show.legend = FALSE  colour="black", 
  # scale_x_discrete(labels=c("Studied","Lure")) +
  scale_fill_manual(values = c("#D2E6BD", "#4AA7B4"),# c("#6798D4", "#D87B49", "#4FA77C"), 
                    labels = c("Violation", "Non-violation")) +
  # scale_fill_manual(values = c("#feb24c", "#91bfdb"), #6798D4 #D87B49 #4FA77C
  # labels = c("Violation", "Non-violation")) +
  coord_cartesian(ylim = c(0, 120), clip = "on") +
  scale_y_continuous(breaks = c(0,25,50,75,100)) +
  labs(x = "Subject Number", y = "Accuracy (%)", fill ="Violation Condition") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position="top")
# legend.position=c(0.8, 0.90))
t2.plot1.sn

2.2.2 Accuracy Analysis: Signal Detection Measure


t2sdtL.1.1 <- t2 %>% filter(cType=="old") %>% mutate( ) %>% group_by(SN, cCong, p_resp, r_resp) %>%
  summarise(n = sum(Corr==1 | Corr==0)/60) %>%
  unite(p_resp.r_resp, c(p_resp,r_resp), sep = "_", remove = TRUE) %>% ungroup() %>%
  spread(p_resp.r_resp, n)
## `summarise()` has grouped output by 'SN', 'cCong', 'p_resp'. You can override using the `.groups` argument.
t2sdtL.1.1 <- plyr::rename(t2sdtL.1.1, c("p.old_r.new" = "Miss", "p.old_r.old" = "Hit"))
t2sdtL.1.1 <- t2sdtL.1.1[,c(1, 2, 4, 3)]
# t2sdtL.1.1 %>% kable(digits=2)

t2sdtL.1.2 <- t2 %>% filter(cType=="new") %>% mutate( ) %>% group_by(SN, cCong, p_resp, r_resp) %>%
  # summarise(n = sum(Corr==1 | Corr==0)/60) %>%
  summarise(n = sum(Corr==1 | Corr==0)/120) %>%
  unite(p_resp.r_resp, c(p_resp,r_resp), sep = "_", remove = TRUE) %>% ungroup() %>%
  spread(p_resp.r_resp, n)
## `summarise()` has grouped output by 'SN', 'cCong', 'p_resp'. You can override using the `.groups` argument.
t2sdtL.1.2 <- plyr::rename(t2sdtL.1.2, c("p.new_r.new" = "CR", "p.new_r.old" = "FA"))
t2sdtL.1.2 <- t2sdtL.1.2[,c(1, 2, 4, 3)]
# t2sdtL.1.2 %>% kable(digits=2)


t2sdtL.3 <- merge(t2sdtL.1.1, t2sdtL.1.2, by="SN")
t2sdtL.3$cCong <- t2sdtL.3$cCong.x
t2sdtL.3[is.na(t2sdtL.3)] <-0
t2sdtL.3[,3] <- t2sdtL.3[,3] + 0.05
t2sdtL.3[,6] <- t2sdtL.3[,6] + 0.05
t2sdtL.3 <- t2sdtL.3 %>% dplyr::select(SN, cCong, Hit, Miss, FA, CR) %>%
  group_by(SN, cCong) %>% dplyr::mutate(mPerf = (Hit - FA)) %>% ungroup()
# t2sdtL.3 %>% kable(digits=2)

indices <- psycho::dprime(t2sdtL.3$Hit, t2sdtL.3$FA, n_targets=1, n_distractors = 1, adjusted=F)

t2sdtL.4 <- cbind(t2sdtL.3, indices)
# t2sdtL.4 %>% kable(digits=2)

t2sdtL.G <- t2sdtL.4 %>%
  dplyr::select(SN, cCong, mPerf, Hit, FA, dprime, aprime, beta, c) %>%
  group_by(cCong) %>%
  summarise(m.p = mean(mPerf), m.d = mean(dprime), m.a = mean(aprime), m.b = mean(beta), m.c = mean(c)) %>% ungroup()
t2sdtL.G %>% kable(digits=2)
cCong m.p m.d m.a m.b m.c
v 0.51 1.62 0.85 2.34 0.50
nv 0.47 1.51 0.84 2.37 0.55


# hit - fa
t2sdtW.1 <- t2sdtL.4 %>% dplyr::select(SN, cCong, mPerf) %>% spread(key= cCong, value=mPerf)
# t2sdtW.1 %>% kable(digits=2)

t2sdtG.1 <- Rmisc::summarySEwithin(data = t2sdtL.4, measurevar = "mPerf", idvar = "SN", withinvars = c("cCong"))
# t2sdtG.1 %>% kable(digits=2)

t2.plot2 <- ggplot(data=t2sdtL.4, aes(x=cCong, y=mPerf, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = 1,
               width = 0.8, show.legend = FALSE) + # , colour="white", size = 1, show.legend = FALSE, colour="black"
  geom_pointrange(data=t2sdtG.1, aes(x = cCong, y=mPerf, ymin = mPerf-se, ymax = mPerf+se),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) + #"darkblue"
  geom_segment(data=t2sdtW.1, inherit.aes = FALSE,
               aes(x=1, y=v,
                   xend=2, yend=nv),
               color="gray60") +
  geom_point(data=t2sdtL.4, position=position_dodge(0.80), color="gray60", size=2, show.legend = FALSE) +
  scale_x_discrete(labels=c("Violation", "Non-violation")) +
  scale_fill_manual(values =  c("#feb24c", "#91bfdb")) +
  coord_cartesian(ylim = c(0, 1), clip = "on") +
  # scale_y_continuous(breaks = c(0,0.25,0.50,75,100)) +
  labs(x = "Violation Condition", y = "Hit - FA") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"),
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))

# d prime
t2sdtW.2 <- t2sdtL.4 %>% dplyr::select(SN, cCong, dprime) %>% spread(key= cCong, value=dprime)
# t2sdtW.2 %>% kable(digits=2)

t2sdtG.2 <- Rmisc::summarySEwithin(data = t2sdtL.4, measurevar = "dprime", idvar = "SN", withinvars = c("cCong"))
# t2sdtG.2 %>% kable(digits=2)

t2.plot3 <- ggplot(data=t2sdtL.4, aes(x=cCong, y=dprime, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = 1,
               width = 0.8, show.legend = FALSE) + # , colour="white", size = 1, show.legend = FALSE, colour="black"
  geom_pointrange(data=t2sdtG.2, aes(x = cCong, y=dprime, ymin = dprime-se, ymax = dprime+se),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) + #"darkblue"
  geom_segment(data=t2sdtW.2, inherit.aes = FALSE,
               aes(x=1, y=v,
                   xend=2, yend=nv),
               color="gray60") +
  geom_point(data=t2sdtL.4, position=position_dodge(0.80), color="gray60", size=2, show.legend = FALSE) +
  scale_x_discrete(labels=c("Violation", "Non-violation")) +
  scale_fill_manual(values =  c("#feb24c", "#91bfdb")) +
  coord_cartesian(ylim = c(0, 4), clip = "on") +
  # scale_y_continuous(breaks = c(0,0.25,0.50,75,100)) +
  labs(x = "Violation Condition", y = "d prime") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"),
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))
# t2.plot3
t2.plot2.all <-ggarrange(t2.plot2, t2.plot3, ncol = 2)
t2.plot2.all

# criterion
t2sdtW.3 <- t2sdtL.4 %>% dplyr::select(SN, cCong, beta) %>% spread(key= cCong, value=beta)
# t2sdtW.3 %>% kable(digits=2)

t2sdtG.3 <- Rmisc::summarySEwithin(data = t2sdtL.4, measurevar = "beta", idvar = "SN", withinvars = c("cCong"))
# t2sdtG.3 %>% kable(digits=2)

t2.plot4 <- ggplot(data=t2sdtL.4, aes(x=cCong, y=beta, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = 1,
               width = 0.8, show.legend = FALSE) + # , colour="white", size = 1, show.legend = FALSE, colour="black"
  geom_pointrange(data=t2sdtG.3, aes(x = cCong, y=beta, ymin = beta-se, ymax = beta+se),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) + #"darkblue"
  geom_segment(data=t2sdtW.3, inherit.aes = FALSE,
               aes(x=1, y=v,
                   xend=2, yend=nv),
               color="gray60") +
  geom_point(data=t2sdtL.4, position=position_dodge(0.80), color="gray60", size=2, show.legend = FALSE) +
  scale_x_discrete(labels=c("Violation", "Non-violation")) +
  scale_fill_manual(values =  c("#feb24c", "#91bfdb")) +
  coord_cartesian(ylim = c(0, 4), clip = "on") +
  # scale_y_continuous(breaks = c(0,0.25,0.50,75,100)) +
  labs(x = "Violation Condition", y = "Beta (response bias)") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"),
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))
t2.plot4


t2accL.1 <- t2sdtL.1.1 # %>% filter(SN!=4) 
t2.dpr.ttest <- t.test(Hit ~ cCong, data=t2accL.1,
                      alternative = c("two.sided"),
                      paired=T, var.equal=F,
                      conf.level=.95)
t2.dpr.ttest
## 
##  Paired t-test
## 
## data:  Hit by cCong
## t = 2.7136, df = 31, p-value = 0.01077
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.009315413 0.065684587
## sample estimates:
## mean of the differences 
##                  0.0375


RM-ANOVA 1


t2accL.1 <- t2sdtL.1.1
t2.acc.aov1 <- aov_ez(id="SN", dv = "Hit", data = t2accL.1, within = c("cCong"))
# test_sphericity(t2.acc.aov1)
# summary(t2.acc.aov1)
# anova(t2.acc.aov1, es = "pes")
nice(t2.acc.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 0.00 7.36 * .192 .011


RM ANOVA 1 - post-hoc


t2.acc.m1 <- emmeans(t2.acc.aov1, pairwise ~ cCong) # adjust="bon" , "tukey"
t2.acc.m1$contrasts %>% kable(digits=3)
contrast estimate SE df t.ratio p.value
v - nv 0.038 0.014 31 2.714 0.011
plot(t2.acc.m1, horizontal = FALSE, comparisons = T)


2.2.3 RT Analysis

tt2 <- t2 %>% filter(Corr == 1, cType=="old") 
tt3 <- tt2 %>% filter(RT > 200 & RT < 10000) %>%
  group_by(SN) %>% # grouping by participants
  nest() %>%
  mutate(lbound = map(data, ~mean(.$RT)-3*sd(.$RT)),
         ubound = map(data, ~mean(.$RT)+3*sd(.$RT))) %>% # make new data (3sd cut)
  unnest(c(lbound, ubound))%>%
  unnest(data) %>%
  mutate(Outlier = (RT < lbound)|(RT > ubound)) %>% # set outlier
  filter(Outlier == FALSE) %>% # filtering outlier
  ungroup()

100-100*(nrow(tt3)/nrow(tt2))
## [1] 2.111324

# subject-level, long format (SN/Btw/Block)
t2rtL <- tt2 %>% group_by(SN, cCong) %>%
  dplyr::summarise(RT = mean(RT)) %>%
  ungroup()
## `summarise()` has grouped output by 'SN'. You can override using the `.groups` argument.
# t2rtL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t2rtW1 <- t2rtL %>% spread(key=cCong, value = RT)
# t2rtW1 %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t2rtG <- t2rtL %>% group_by(cCong) %>%
  summarise(RT.m = mean(RT), RT.sd = sd(RT)) %>%
  ungroup()
t2rtG$RT.se <- Rmisc::summarySEwithin(data = t2rtL, measurevar = "RT", 
                                        idvar = "SN", withinvars = c("cCong"))$se
t2rtG$RT.ci <- Rmisc::summarySEwithin(data = t2rtL, measurevar = "RT", 
                                        idvar = "SN", withinvars = c("cCong"))$ci
t2rtG <- t2rtG %>% 
  mutate(lower.ci = RT.m-RT.ci,
         upper.ci = RT.m+RT.ci)
# for between-subject design. check help(summarySE)
t2rtG %>% kable(digits=2)
cCong RT.m RT.sd RT.se RT.ci lower.ci upper.ci
v 1573.0 397.24 27.93 56.97 1516.03 1629.98
nv 1591.6 403.22 27.93 56.97 1534.63 1648.58


t2.rt.plot1 <- ggplot(data=t2rtL, aes(x=cCong, y=RT, fill=cCong)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02, show.legend = FALSE ) + # , show.legend = FALSE  colour="black", 
  geom_point(data=t2rtL, aes(x=cCong, y=RT, fill=cCong), 
             position=position_dodge(width = 0.8), 
             size=2, show.legend = FALSE, color="gray60") + # , 
  geom_segment(data=t2rtW1, inherit.aes = FALSE,
               aes(x=1, y=v,
                   xend=2, yend=nv),
               color="gray60") +
  geom_pointrange(data=t2rtG, aes(x = cCong, y=RT.m, ymin = lower.ci, ymax=upper.ci),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) +
  scale_x_discrete(labels=c("Violation", "Non-violation")) +
  scale_fill_manual(values = c("#D2E6BD", "#4AA7B4")) +# c("#6798D4", "#D87B49", "#4FA77C"), 
  # scale_fill_manual(values = c("#feb24c", "#91bfdb"), #6798D4 #D87B49 #4FA77C
  # labels = c("Violation", "Non-violation")) +
  coord_cartesian(ylim = c(0, 3000), clip = "on") +
  # scale_y_continuous(breaks = c(0,25,50,75,100)) +
  labs(x = "Stimuli Type", y = "Response Time (ms)", fill ="Violation Condition") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position="top")
# legend.position=c(0.8, 0.90))
t2.rt.plot1


Paired T-Test


t2.rt.ttest <- t.test(RT ~ cCong, data=t2rtL,
                       alternative = c("two.sided"),
                       paired=T, var.equal=F,
                       conf.level=.95)
t2.rt.ttest
## 
##  Paired t-test
## 
## data:  RT by cCong
## t = -0.40773, df = 31, p-value = 0.6863
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -111.63683   74.43789
## sample estimates:
## mean of the differences 
##               -18.59947




2.3 Relating priming effect to memory performance


# categorization
t1 <- read.csv("data/reactMEM_p1_all.csv", header = T)
t1$SN = factor(t1$SN)
t1$Block = factor(t1$Block)
t1$cCong = factor(t1$cCong, levels=c(1,2), labels=c("v","nv")) 
t1$cnRep = factor(t1$cnRep, levels=c(1,2,3,5,4), labels=c(1,2,3,4,0))
t1$cnPair = factor(t1$cnPair)
t1$cIdt = factor(t1$cIdt)
t1$cCat = factor(t1$cCat)
t1$Im_Name = factor(t1$Im_Name)
t1$Im_Idx = factor(t1$Im_Idx)
t1$RT <- t1$RT*1000;
t1$Corr <- as.numeric(t1$Corr==1)

# recognition
t2 <- read.csv("data/reactMEM_p2_all.csv", header = T)
t2$SN = factor(t2$SN)
t2$cType = factor(t2$cType, levels=c(1,2), labels=c("old","new"))
t2$cCong = factor(t2$cCong, levels=c(1,2,999), labels=c("v","nv","new")) 
t2$cIdt = factor(t2$cIdt)
t2$dIm_Name = factor(t2$dIm_Name)
t2$Resp <- factor(t2$Resp);
t2$mResp <- factor(t2$Resp, levels=c(1,2,3,4), labels=c("new", "new", "old", "old"))
t2$RT <- t2$RT*1000;
t2$Corr <- as.numeric(t2$Corr==1)
t2$Conf <- as.double(t2$Conf)
t2$p_resp = factor(t2$cType, levels=c('old','new'), labels=c("p.old","p.new"))
t2$r_resp = factor(t2$mResp, levels=c('old','new'), labels=c("r.old", "r.new"))

# RT preprocessing
t3.1 <- t1 %>% filter(Corr == 1) %>% # filter(RT > 200 & RT < 10000) %>%
  group_by(SN) %>% # grouping by participants
  nest() %>%
  mutate(lbound = map(data, ~mean(.$RT)-3*sd(.$RT)),
         ubound = map(data, ~mean(.$RT)+3*sd(.$RT))) %>% # make new data (3sd cut)
  unnest(c(lbound, ubound))%>%
  unnest(data) %>%
  mutate(Outlier = (RT < lbound)|(RT > ubound)) %>% # set outlier
  filter(Outlier == FALSE) %>% # filtering outlier
  ungroup()
t3.1 <- t3.1 %>% filter(cnRep != 1, cnRep != 2)
t3.2 <- t3.1 %>% filter(cnPair == 2) %>% 
  select(SN, cCong, cnRep, cIdt, Im_Idx, RT)
t3.3 <- t3.2 %>% group_by(SN) %>% 
  spread(key = cnRep, value = RT) %>% ungroup()
t3.4 <- t3.3 %>% group_by(SN, cIdt) %>% mutate(prime = `4` - `3`)

# Acc preprocessing
t4.1 <- t2 %>% filter(cType == 'old') %>% 
  select (SN, cCong, cIdt, dIm_Name, Corr, RT, Conf)

t5 <- merge(t4.1, t3.4, key=c(SN,cCong,cIdt))


2.3.1 Logistic Regression


p1.glm <- t5 #  %>% filter(SN!=4)
p1.glm <- na.omit(p1.glm)
table(p1.glm$SN,p1.glm$cCong)
##     
##       v nv new
##   1  46 57   0
##   2  53 53   0
##   3  55 54   0
##   4  50 53   0
##   5  56 59   0
##   6  50 56   0
##   7  45 51   0
##   8  56 54   0
##   9  53 55   0
##   10 46 46   0
##   11 49 53   0
##   12 33 39   0
##   13 54 54   0
##   14 46 47   0
##   15 55 53   0
##   16 51 53   0
##   17 53 53   0
##   18 46 52   0
##   19 55 57   0
##   20 58 54   0
##   21 57 58   0
##   22 55 55   0
##   23 53 52   0
##   24 48 52   0
##   25 54 58   0
##   26 43 48   0
##   27 56 55   0
##   28 57 52   0
##   29 56 52   0
##   30 52 52   0
##   31 50 54   0
##   32 58 56   0
table(p1.glm$Corr,p1.glm$cCong)
##    
##       v  nv new
##   0 724 811   0
##   1 925 886   0
glimpse(p1.glm, width=70)
## Rows: 3,346
## Columns: 12
## $ SN       <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
## $ cCong    <fct> nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv,…
## $ cIdt     <fct> 10, 100, 101, 106, 107, 108, 114, 115, 119, 120, 17…
## $ dIm_Name <fct> 144_2_out.jpg, 136_2_out.jpg, 155_2_out.jpg, 101_2_…
## $ Corr     <dbl> 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, …
## $ RT       <dbl> 1211.4, 969.5, 982.7, 2727.3, 1109.6, 1099.3, 1000.…
## $ Conf     <dbl> 1, 1, 2, 1, 2, 2, 2, 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, …
## $ Im_Idx   <fct> 144, 136, 155, 101, 194, 132, 3, 172, 32, 14, 22, 1…
## $ `3`      <dbl> 529.5, 568.3, 538.0, 479.2, 536.1, 506.8, 478.5, 51…
## $ `4`      <dbl> 479.1, 531.8, 670.9, 494.6, 526.4, 467.2, 457.5, 43…
## $ `0`      <dbl> 481.6, 538.8, 470.1, 521.0, 538.9, 536.4, 587.7, 44…
## $ prime    <dbl> -50.4, -36.5, 132.9, 15.4, -9.7, -39.6, -21.0, -80.…


p1.glm.v <- p1.glm %>% filter(cCong == "v")
p1.glm1 <- glm(Corr ~ prime, data=p1.glm.v, family=binomial()) 
# summary(p1.glm1)  %>% kable(digits=2)
# sjPlot::tab_model(p1.glm1)

p1.glm1.2 <- glm(Corr ~ `0`, data=p1.glm.v, family=binomial())
# summary(p1.glm1.2) # display results
# anova(p1.glm1.2)
# sjPlot::tab_model(p1.glm1.2)

p1.glm.v$prime2 <- p1.glm.v$prime^2
p1.glm1.1 <- glm(Corr ~ prime + prime2, data=p1.glm.v, family=binomial()) 
# summary(p1.glm1.1) # display results
# anova(p1.glm1.1)
# sjPlot::tab_model(p1.glm1.1)
sjPlot::tab_model(p1.glm1)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.27 1.16 – 1.41 <0.001
prime 1.00 1.00 – 1.00 0.240
Observations 1649
R2 Tjur 0.001
sjPlot::tab_model(p1.glm1.1)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.30 1.17 – 1.44 <0.001
prime 1.00 1.00 – 1.00 0.218
prime2 1.00 1.00 – 1.00 0.258
Observations 1649
R2 Tjur 0.002
sjPlot::tab_model(p1.glm1.2)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.63 1.20 – 2.22 0.002
0 1.00 1.00 – 1.00 0.099
Observations 1649
R2 Tjur 0.002


p3_logreg.v <- ggplot(p1.glm.v, aes(x=prime, y=Corr)) + 
  geom_point(alpha=.5, size=2) +
  # ggbeeswarm::geom_quasirandom(dodge.width = 1, color = "blue", size = 1, alpha = 0.2,
                               # show.legend = FALSE) +
  stat_smooth(method="glm", se=T, method.args = list(family=binomial)) + 
  labs(x = "Diff Score (final B rt - 3rd B rt)", y = "Correctness") +
  scale_y_continuous(breaks = c(0, 0.5, 1)) +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))
# p3_logreg.v

p3_logreg.v.2 <- ggplot(p1.glm.v, aes(x=`0`, y=Corr)) + 
  geom_point(alpha=.5, size=2) +
  # ggbeeswarm::geom_quasirandom(dodge.width = 1, color = "blue", size = 1, alpha = 0.2,
  # show.legend = FALSE) +
  stat_smooth(method="glm", se=T, method.args = list(family=binomial)) + 
  labs(x = "B` RT", y = "Correctness") +
  scale_y_continuous(breaks = c(0, 0.5, 1)) +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))


p3_logreg.plot <- ggarrange(p3_logreg.v, p3_logreg.v.2, ncol = 2,
          labels=c("A) Violation, Prime", "B) Violation, B`"))
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
p3_logreg.plot


p1.glm.nv <- p1.glm %>% filter(cCong == "nv")
p1.glm.nv <- p1.glm.nv %>% group_by(SN) %>% mutate(prime_z = (prime-mean(prime))/sd(prime)) %>% ungroup()

p1.glm2 <- glm(Corr ~ prime, data=p1.glm.nv, family=binomial()) 
# summary(p1.glm2) # display results
# anova(p1.glm2)

p1.glm.nv$prime2 <- p1.glm.nv$prime^2
p1.glm2.1 <- glm(Corr ~ prime + prime2, data=p1.glm.nv, family=binomial()) 
# summary(p1.glm2.1) # display results
# anova(p1.glm2.1)

p1.glm2.2 <- glm(Corr ~ `0`, data=p1.glm.nv, family=binomial())
# summary(p1.glm2.2) # display results
# anova(p1.glm2.2)

sjPlot::tab_model(p1.glm2)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.09 0.99 – 1.20 0.067
prime 1.00 1.00 – 1.00 0.753
Observations 1697
R2 Tjur 0.000
sjPlot::tab_model(p1.glm2.1)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.09 0.98 – 1.21 0.096
prime 1.00 1.00 – 1.00 0.764
prime2 1.00 1.00 – 1.00 0.872
Observations 1697
R2 Tjur 0.000
sjPlot::tab_model(p1.glm2.2)
  Corr
Predictors Odds Ratios CI p
(Intercept) 1.61 1.20 – 2.17 0.002
0 1.00 1.00 – 1.00 0.007
Observations 1697
R2 Tjur 0.004
# hist(residuals(p1.glm1, type="deviance"))
# qqnorm(residuals(p1.glm1, type="deviance"))


# plotting
newdat.nv <- data.frame(prime=seq(min(p1.glm.nv$prime), max(p1.glm.nv$prime),len=dim(p1.glm.nv)[1]))
newdat.nv$Corr = predict(p1.glm2, newdata=p1.glm.nv, type="response")

p3_logreg.nv <- ggplot(p1.glm.nv, aes(x=prime, y=Corr)) + 
  geom_point(alpha=.5, size=2) +
  # ggbeeswarm::geom_quasirandom(dodge.width = 1, color = "blue", size = 1, alpha = 0.2,
  # show.legend = FALSE) +
  stat_smooth(method=glm, se=T, method.args = list(family=binomial)) + 
  labs(x = "Diff Score (final B rt - 3rd B rt)", y = "Correctness") +
  scale_y_continuous(breaks = c(0, 0.5, 1)) +
  # ggtitle("Non-violation Condition") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))
# p3_logreg.nv

p3_logreg.nv.2 <- ggplot(p1.glm.nv, aes(x=`0`, y=Corr)) + 
  geom_point(alpha=.5, size=2) +
  # ggbeeswarm::geom_quasirandom(dodge.width = 1, color = "blue", size = 1, alpha = 0.2,
  # show.legend = FALSE) +
  stat_smooth(method="glm", se=T, method.args = list(family=binomial)) + 
  labs(x = "B` RT", y = "Correctness") +
  scale_y_continuous(breaks = c(0, 0.5, 1)) +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position=c(0.8, 0.85))
# p3_logreg.nv.2

p3_logreg.plot <- ggarrange(p3_logreg.nv, p3_logreg.nv.2, ncol = 2,
          labels=c("A) non-Violation, Prime", "B) non-Violation, B`"))
## `geom_smooth()` using formula 'y ~ x'
## `geom_smooth()` using formula 'y ~ x'
p3_logreg.plot


2.3.2 Generalized Linear Mixed Effect Model


p1.glm.v <- p1.glm.v %>% group_by(SN) %>% 
  mutate(prime_z = (prime-mean(prime))/sd(prime)) %>% ungroup()

(nc <- detectCores())
## [1] 8
cl <- makeCluster(rep("localhost", nc))
p1.glm <- p1.glm %>% group_by(SN) %>%
  mutate(prime_z_a = (prime-mean(prime))/sd(prime)) %>% ungroup()
p1.glm.mixed <- afex::mixed(Corr ~  cCong*prime + (1|SN) + (1|Im_Idx),
                              data = p1.glm,
                              family = binomial(link="logit"),
                              method = "LRT", cl = cl,
                              expand_re = T,
                              control = lmerControl(optimizer = "bobyqa",
                                                    optCtrl = list(maxfun = 1e7)))
## Fitting 4 (g)lmer() models.
p1.glm.mixed.v <- afex::mixed(Corr ~ prime_z + (1|SN) + (1|Im_Idx),
                              data = p1.glm.v,
                              family = binomial(link="logit"),
                              method = "LRT", cl = cl,
                              expand_re = T,
                              control = lmerControl(optimizer = "bobyqa",
                                                    optCtrl = list(maxfun = 1e7)))
## Fitting 2 (g)lmer() models.
p1.glm.mixed.nv <- afex::mixed(Corr ~ prime_z + (1|SN) + (1|Im_Idx),
                               data = p1.glm.nv,
                               family = binomial(link="logit"),
                               method = "LRT", cl = cl,
                               expand_re = T,
                               control = lmerControl(optimizer = "bobyqa",
                                                     optCtrl = list(maxfun = 1e7)))
## Fitting 2 (g)lmer() models.
stopCluster(cl)


sjPlot::tab_model(p1.glm.mixed)
  Corr
Predictors Estimates CI p
(Intercept) 0.28 0.00 – 0.56 0.046
cCong [nv] -0.17 -0.32 – -0.02 0.031
prime 0.00 -0.00 – 0.00 0.435
cCong [nv] * prime -0.00 -0.00 – 0.00 0.740
Random Effects
σ2 3.29
τ00 Im_Idx 0.30
τ00 SN 0.50
ICC 0.20
N SN 32
N Im_Idx 240
Marginal R2 / Conditional R2 0.002 / 0.197
sjPlot::tab_model(p1.glm.mixed.v)
  Corr
Predictors Estimates CI p
(Intercept) 0.31 0.04 – 0.58 0.027
prime_z 0.05 -0.06 – 0.16 0.387
Random Effects
σ2 3.29
τ00 Im_Idx 0.36
τ00 SN 0.47
ICC 0.20
N SN 32
N Im_Idx 240
Marginal R2 / Conditional R2 0.001 / 0.202
sjPlot::tab_model(p1.glm.mixed.nv)
  Corr
Predictors Estimates CI p
(Intercept) 0.11 -0.15 – 0.38 0.409
prime_z 0.06 -0.05 – 0.16 0.271
Random Effects
σ2 3.29
τ00 Im_Idx 0.17
τ00 SN 0.48
ICC 0.17
N SN 32
N Im_Idx 240
Marginal R2 / Conditional R2 0.001 / 0.166


2.3.3 Binning and ANOVA


p1.bin <- p1.glm  %>%  group_by(SN, cCong) %>% # grouping by participants
  mutate(med = median(prime)) %>% # make new data (3sd cut)
  mutate(prime_bin = (prime >= med)) %>% # set outlier
  ungroup()

p1.bin$prime_bin <- as.numeric(p1.bin$prime_bin == T)
p1.bin$prime_bin <- p1.bin$prime_bin+1
p1.bin$prime_bin <- factor(p1.bin$prime_bin, levels=c(1,2), labels = c("low", "high"))
glimpse(p1.bin)
## Rows: 3,346
## Columns: 15
## $ SN        <fct> 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, …
## $ cCong     <fct> nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, nv, …
## $ cIdt      <fct> 10, 100, 101, 106, 107, 108, 114, 115, 119, 120, 17, 18, 19,…
## $ dIm_Name  <fct> 144_2_out.jpg, 136_2_out.jpg, 155_2_out.jpg, 101_2_in.jpg, 1…
## $ Corr      <dbl> 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, 1, 0, 1, 1, …
## $ RT        <dbl> 1211.4, 969.5, 982.7, 2727.3, 1109.6, 1099.3, 1000.0, 799.1,…
## $ Conf      <dbl> 1, 1, 2, 1, 2, 2, 2, 2, 2, 1, 1, 2, 1, 2, 1, 2, 1, 1, 1, 2, …
## $ Im_Idx    <fct> 144, 136, 155, 101, 194, 132, 3, 172, 32, 14, 22, 141, 5, 15…
## $ `3`       <dbl> 529.5, 568.3, 538.0, 479.2, 536.1, 506.8, 478.5, 514.4, 538.…
## $ `4`       <dbl> 479.1, 531.8, 670.9, 494.6, 526.4, 467.2, 457.5, 433.8, 509.…
## $ `0`       <dbl> 481.6, 538.8, 470.1, 521.0, 538.9, 536.4, 587.7, 448.5, 482.…
## $ prime     <dbl> -50.4, -36.5, 132.9, 15.4, -9.7, -39.6, -21.0, -80.6, -28.9,…
## $ prime_z_a <dbl> -0.55233883, -0.40079548, 1.44607094, 0.16503904, -0.1086111…
## $ med       <dbl> 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, …
## $ prime_bin <fct> low, low, high, high, low, low, low, low, low, high, low, lo…
# subject-level, long format (SN/Btw/Block)
t5rtL <- p1.bin %>% group_by(SN, cCong, prime_bin) %>%
  dplyr::summarise(p2Acc = mean(Corr)*100) %>%
  ungroup()
# t5rtL %>% kable(digits=2)

# subject-level, wide format (SN/Btw/Block)
t5rtW1 <- t5rtL %>% spread(key=prime_bin, value = p2Acc)
# t5rtW1 %>% kable(digits=2)

# summary table: grand mean (eyerep/locrep)
t5rtG <- t5rtL %>% group_by(cCong, prime_bin) %>%
  summarise(p2Acc.m = mean(p2Acc), p2Acc.sd = sd(p2Acc)) %>%
  ungroup()
t5rtG$p2Acc.se <- Rmisc::summarySEwithin(data = t5rtL, measurevar = "p2Acc", 
                                      idvar = "SN", withinvars = c("cCong","prime_bin"))$se
t5rtG$p2Acc.ci <- Rmisc::summarySEwithin(data = t5rtL, measurevar = "p2Acc", 
                                      idvar = "SN", withinvars = c("cCong","prime_bin"))$ci
t5rtG <- t5rtG %>% 
  mutate(lower.ci = p2Acc.m-p2Acc.ci,
         upper.ci = p2Acc.m+p2Acc.ci)
# for between-subject design. check help(summarySE)
t5rtG %>% kable(digits=2)
cCong prime_bin p2Acc.m p2Acc.sd p2Acc.se p2Acc.ci lower.ci upper.ci
v low 55.71 18.44 1.78 3.64 52.07 59.35
v high 56.78 16.06 1.22 2.50 54.28 59.27
nv low 51.99 18.71 1.56 3.18 48.81 55.17
nv high 52.59 17.92 1.55 3.15 49.43 55.74


t5.plot1 <- ggplot(data=t5rtL, aes(x=cCong, y=p2Acc, fill=prime_bin)) +
  stat_summary(fun = mean, geom = "bar", position="dodge", na.rm = TRUE, alpha = .9, 
               width = 0.8,size = 1.02) + # , show.legend = FALSE  colour="black", 
  geom_point(data=t5rtL, aes(x=cCong, y=p2Acc, fill=prime_bin), 
             position=position_dodge(width = 0.8), 
             size=2, show.legend = FALSE, color="gray60") + # , 
  geom_segment(data=filter(t5rtW1, cCong == "v"), inherit.aes = FALSE,
               aes(x=0.8, y=filter(t5rtW1, cCong == "v")$low,
                   xend=1.2, yend=filter(t5rtW1, cCong == "v")$high),
               color="gray60") +
  
  geom_segment(data=filter(t5rtW1, cCong == "nv"), inherit.aes = FALSE,
               aes(x=1.8, y=filter(t5rtW1, cCong == "nv")$low,
                   xend=2.2, yend=filter(t5rtW1, cCong == "nv")$high),
               color="gray60") +
  
  geom_pointrange(data=t5rtG, aes(x = cCong, y=p2Acc.m, ymin = lower.ci, ymax = upper.ci),
                  position = position_dodge(0.80), color = "darkred", size = 1, show.legend = FALSE) +
  scale_x_discrete(labels=c("Violation","Non-violation")) +
  scale_fill_manual(values = c("#D2E6BD", "#4AA7B4"),# c("#6798D4", "#D87B49", "#4FA77C"), 
                    labels = c("Low", "High")) +
  # scale_fill_manual(values = c("#feb24c", "#91bfdb"), #6798D4 #D87B49 #4FA77C
  # labels = c("Violation", "Non-violation")) +
  coord_cartesian(ylim = c(0, 120), clip = "on") +
  scale_y_continuous(breaks = c(0,25,50,75,100)) +
  labs(x = "Violation Condition", y = "p2 Recognition Hit (%)", fill ="p1 4th-3rd B rt") +
  # ggtitle("Association Memory Accuracy") +
  theme_bw(base_size = 15) +
  theme(axis.title = element_text(face = "bold", size = 16, color = "black"),
        axis.text = element_text(face = "plain", hjust = 0.5, size = 15, color = "black"),
        axis.line=element_line(),
        strip.text.x = element_text(face = "plain", size = 15, color = "black"),
        strip.background = element_blank(),
        panel.grid.major = element_blank(),
        panel.grid.minor = element_blank(),
        panel.border = element_blank(),
        panel.spacing=unit(1, "lines"),
        # aspect.ratio = 0.2,
        plot.margin = margin(0.3, 0.3, 0.3, 0.3, "cm"), 
        # plot.title = element_text(hjust = 0.5),
        legend.position="top")
# legend.position=c(0.8, 0.90))
t5.plot1


RM-ANOVA


t5.acc.aov1 <- aov_ez(id="SN", dv = "p2Acc", data = t5rtL, within = c("cCong", "prime_bin"))
# test_sphericity(t2.acc.aov1)
# summary(t5.acc.aov1)
# anova(t5.acc.aov1, es = "pes")
nice(t5.acc.aov1, es="pes") %>% kable(digits=3)
Effect df MSE F pes p.value
cCong 1, 31 76.33 6.55 * .174 .016
prime_bin 1, 31 84.81 0.26 .008 .614
cCong:prime_bin 1, 31 92.33 0.02 <.001 .890


RM ANOVA - post-hoc


t5.acc.m1 <- emmeans(t5.acc.aov1, pairwise ~ prime_bin | cCong) # adjust="bon" , "tukey"
t5.acc.m1$contrasts %>% kable(digits=3)
contrast cCong estimate SE df t.ratio p.value
low - high v -1.065 2.328 31 -0.457 0.651
low - high nv -0.593 2.378 31 -0.249 0.805
plot(t5.acc.m1, horizontal = FALSE, comparisons = T)


Paired T-Test


Violation


t5rtL.v <- t5rtL %>% filter(cCong=="v")
t5.acc.ttest <- t.test(p2Acc ~ prime_bin, data=t5rtL.v,
                       alternative = c("two.sided"),
                       paired=T, var.equal=F,
                       conf.level=.95)
t5.acc.ttest
## 
##  Paired t-test
## 
## data:  p2Acc by prime_bin
## t = -0.45735, df = 31, p-value = 0.6506
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -5.812061  3.682881
## sample estimates:
## mean of the differences 
##                -1.06459


Non-Violation


t5rtL.nv <- t5rtL %>% filter(cCong=="nv")

t5.acc.ttest <- t.test(p2Acc ~ prime_bin, data=t5rtL.nv,
                       alternative = c("two.sided"),
                       paired=T, var.equal=F,
                       conf.level=.95)
t5.acc.ttest
## 
##  Paired t-test
## 
## data:  p2Acc by prime_bin
## t = -0.24923, df = 31, p-value = 0.8048
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -5.441678  4.256527
## sample estimates:
## mean of the differences 
##              -0.5925753